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ABSTRACT 

A new calculation of the Galactic diffuse gamma-ray spectrum from the decay of 
secondary particles produced by interactions of cosmic-ray protons with interstellar 
matter is presented. The calculation utilizes the modern Monte Carlo event generators, 
Hadrin, Fritiof and Pythia, which simulate high-energy proton-proton collisions 
and are widely used in studies of nuclear and particle physics, in addition to scaling 
calculation. This study is motivated by the result on the Galactic diffuse gamma-ray 
flux observed by the EGRET detector on the Compton Gamma-ray Observatory, 
which indicates an excess above about 1 GeV of the observed intensity compared with 
a model prediction. The prediction is based on cosmic-ray interactions with interstellar 
matter, in which secondary pion productions are treated by a simple model. With 
the improved interaction model used here, however, the diffuse gamma-ray flux agrees 
rather well with previous calculations within uncertainties, which mainly come from 
the unobservable demodulated cosmic-ray spectrum in interstellar space. As a possible 
solution to the excess flux, flatter spectra of cosmic-ray protons have been tested and 
we found that the power-law spectrum with an index of about —(2.4 ~ 2.5) gives a 
better fit to the EGRET data, though the spectrum is not explained completely. 



Subject headings: cosmic rays — gamma rays — interstellar medium 



1. Introduction 

Observations of diffuse gamma-ray emission from the Galactic plane give us some knowledge 
on Galactic cosmic-rays, interstellar medium and interaction between them. The observed features, 
most recently by the Energetic Gamma Ray Experiment Telescope (EGRET) on the Compton 
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Gamma-Ray Observatory, are described fairly well by a model based on dynamic balance and 
realistic interstellar matter and photon distributions (Hunter et al. 1996). 

However, the observed intensity exceeds the model prediction by as much as 60% for energies 
above about 1 GeV. One of the possible explanation of this discrepancy is that the theory of 
diffuse gamma-ray production may not be adequate at high energies: astrophysical gamma-rays 
above 1 GeV have never been measured before with high statistical accuracy (Hunter et al. 1996). 

There are three basic components for the production of Galactic diffuse gamma-rays: nuclear 
interactions between cosmic-rays and matter, bremsstrahlung collisions between electrons and 
matter, and inverse Compton scattering of electrons with low-energy photons. Above about 200 
MeV, the first component, more specifically the gamma-ray production from the decay of neutral 
pions produced in cosmic-ray (protons and nuclei) collisions with interstellar matter, is known to 
be the dominant one. Previous works on this component used isobaric models and scaling models 
for the nuclear interaction: see Dermer (1986a) for detail. These studies, however, were focused on 
lower energies. This is natural since the gamma-ray spectrum from pion decay has a peak around 
70 MeV and drops rapidly toward higher energies. In order to compare the observed high-energy 
diffuse gamma-ray emission with model predictions, however, it is necessary to use more detailed 
models which describe high-energy proton-proton (p-p) collisions more accurately. Here we use 
Monte Carlo event generators that are commonly used in the analysis of high-energy physics 
experiments to simulate such collisions. We apply these results for the calculation of the Galactic 
diffuse gamma-ray flux from cosmic-ray interactions up to 10 GeV and compare it with previous 
studies and data from the EGRET detector. 



2. Gamma-ray production in cosmic-ray interaction 

The gamma-ray production source function from nuclear interaction, or the gamma-ray 
spectrum resulting from pion decay in p-Hydrogen atom collisions, is given by Stecker (1970) 

q{E y ) = Aim H dE n ■ 2 / dT p == — 

JE-.+ml/AE-, JTf(T v ) y y/E%-ml dT w 

roo 

= Atthh / dT p j p (Tp)((a n (T p )} 

where the cosmic-ray proton flux is denoted by j p (T p ) [cm -2 s -1 sr -1 GeV -1 ], da(T p , T 7r )/dT 7r is 
the differential cross section for the production of a tt° with kinetic energy T w in the Galactic 
(rest) system due to a collision of a cosmic-ray proton of kinetic energy T p with a H atom at 
rest, E w is the total pion energy and m w is its mass, and nu is the atomic hydrogen density. 
T™ in (T 7r ) is the minimum proton kinetic energy that contributes to the production of a pion with 
energy T w . It can be calculated from the kinematics easily. Following Dermer (1986a), we write 
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da(T p ,T n )/dT w = ((a n (T p ))dN(T p ,T n )/dT n , where ((^(Tp)) is the inclusive cross section for the 
reaction p + p — > 7r° + anything (£ is the 7r° multiplicity) and dA^ (T p , T n )/dT n is the normalized 7r° 
spectrum: 



We studied the three Monte Carlo programs Hadrin (Hanssgen and Ranft 1986), Pythia 
(Stostrand 1994), and Fritiof (Pi 1992), in comparison with the scaling model of the 7r° 
production cross section by p-p interaction given in an analytical form by Stephens and Badhwar 
(1981). This scaling model was used in some of the previous calculations (Stephens and Badhwar 
1981; Dermer 1986a). 

The first program, Hadrin, focuses on a threshold and resonance behavior of inelastic 
hadron-nucleon interactions and uses tabulated cross sections for many possible reaction channels 
based on experimental data. It focuses on describing the nuclear collision at laboratory energies 
below 5 GeV. Fig. [I] shows the comparison of the simulated 7r° kinetic energy spectrum with 
experimental data and analytic calculations based on the isobar model (Stecker 1970) and the 
scaling model (see also Dermer 1986a) in the center-of-mass system (CMS). Also shown in Fig. [2| 
is the same data but in laboratory system (LS). One can see both Hadrin and the isobar model 
calculations show good agreement with experimental data within limited statistical accuracy of 
the data. 

The program called Pythia^] focuses on high-energy p-p colliders and its performance is 
reported in detail (Sjostrand and van Zijl 1987). This is based on a string- fragmentation model, 
but incorporates low transverse momentum (pr) interactions, which is our main concern. Its 
usefulness for fixed target experiments, the situation similar to ours, is not fully studied, however. 
The other program, FRITIOF,^ implements a model for low-p^ hadron-hadron, hadron-nucleus 
and nucleus-nucleus reactions which treats a hadron as a string-like object with its color force 
field stretching like a vortex line. A detailed description of the model and its comparison with 
experimental data are reported (Andersson, Gustafson and Pi 1993). 

Here we compare their predictions with the p-p collision data: due to the experimental 
difficulty in obtaining tt° spectra in high- multiplicities, we use gamma-ray spectra. 

Fig. H shows the CMS rapidity (y* = | In j^rzpS where E* is the pion energy in CMS and 
p* L the longitudinal momentum in CMS) distributions for beam momenta of 12.4, 205, and 300 
GeV/c. Data are obtained in bubble chamber experiments at Argonne (Jaeger et al. 1975a) 



3 We used Pythia version 5.718 with Jetset version 7.408. 

4 We used Fritiof version 7.02 with Ariadne version 4.02R and Jetset version 7.3. 
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2.1. 



Secondary pion production 
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and Fermilab (Jaeger et al. 1975b; Sheng et al. 1975). Predictions and data are normalized to 
match the total cross sections given by the fitted formula given by Dermer (1986b). The scaling 
calculation has been executed as a combination of the analytic Stephens and Badhwar (1981) 
model and a Monte Carlo program Decay (Hanssgen and Ritter 1984) to treat tt° decays. One 
sees a general agreement of data with Monte Carlo predictions, but some overprediction is seen in 
the forward and backward (\y*\ ^ 2) directions at high energies. This tendency is more evident 
for Pythia: it may partly due to its main intention to describe high transverse momentum (pr) 
phenomena, while Fritiof is a model to describe low px hadronic reactions. We do not use 
Pythia to calculate the gamma-ray flux since the forward region plays an importrant role at high 
energies (see below). 

Fig. U shows the comparison of the model prediction of the invariant cross section of tt° 
production with experimental data at ISR. The differential cross sections of gamma-rays based on 
Monte Carlo programs have been converted to invariant cross sections of 7r°'s using Sternheimer 
relation (Sternheimer 1955) following the strategy given in Stephens and Badhwar (1981). Q Monte 
Carlo calculations agree well with the scaling calculations given by Stephens and Badhwar, but 
both of them overpredict the cross sections at higher rapidities, as was seen in the previous figure 
(see the discussion of Stephens and Badhwar about the significance of this discrepancy). This 
results in systematic uncertainties of the final gamma-ray flux at higher energies {E~ ;> 100 GeV). 

2.2. Inclusive cross section 

Measured inclusive cross sections, (Ccr^Tp)), for the reaction p + p — > 7r° + anything are 
well fitted by the analytic parameterization given by Dermer (1986b) for T p < 1 TeV. In order 
to extend the calculation to higher energies, inclusive cross sections have been computed with a 
help of Monte Carlo simulators. In order to match the experimental data, inclusive cross sections 
of 7T° production are computed from inclusive cross sections of gamma-ray production assuming 
all gamma-rays are from 7r° decays. Fig. || compares experimental results with the Dermer's 
parameterization and the Pythia and Fritiof prediction. Total p-p cross sections are calculated 
by a parameterization inspired by Regge theory (Donnachie and Landshoff 1992) in Pythia and 
by a Block-Cahn fit (Block and Cahn 1987) in Fritiof. We adopt Dermer's parameterization 
below 1 TeV/c and use 163(s/1876 GeV 2 ) - 21 mb as a smooth fit to the Pythia prediction above 
1 TeV/c, where s = 2m p {T p + 2m p ) with m p the proton mass, since it is based on newer fit to 
the total cross section and charged particle multiplicity distributions are described well up to 
y/s = 900 GeV (Sjostrand and van Zijl 1987). 



5 We have used the correct formula: £(d 3 a/dp 3 ) = -(l/2)(0/0p)(dV/(dpdfi)) ~ (1/2) (B + 2C)(d 2 a/(dpdfi)) if 
d 2 cr/(dpdfi) oc exp(-Bp - C7p 2 ). 
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2.3. Cosmic-ray proton flux 

Galactic diffuse gamma-rays are mostly produced in the central region of the Galaxy and we 
should use the cosmic-ray spectrum there: however, since our knowledge is limited within the local 
solar neighborhood, the best we can do is to use the local cosmic-ray spectrum which is estimated 
by using appropriate assumptions from the solar-modulated spectrum actually observed. Fig. || 
summarizes the observations and demodulated spectra calculated by different authors. There is 
a factor of about 5 uncertainty in the proton spectrum at T p = 1 GeV, but as we will see these 
low energy protons do not contribute much to the diffuse flux, as has been pointed out before 
(e.g. Stephens and Badhwar 1981). In order to quantify the spectrum uncertainty, we will calculate 
the diffuse gamma-ray flux for three extreme cases: the "maximal" ("minimal") flux corresponds 
to the upper (lower) envelope of calculated demodulated spectra and connects to E~ 2 - 75 spectra 
above 100 GeV, and the "median" flux takes the form of 



JpyTp) 



1.67p p 



, , / 2.5GeV/c 
\ Pp 

6.65 x 10- 6 (i^iv) ~ 2 ' " (E p > 100 GeV) 



-1/2 



(E p < 100 GeV) 

(3) 



where E p = T p + m p , p p = E 2 — m 2 and the unit is cm s sr . These are inspired by 
the parameterization on the demodulated spectrum given by Ormes and Protheroe (1983) for 
E p < 100 GeV and the fit by Honda et al. (1995) for higher energies. The cosmic-ray proton fluxes 
used in this work are shown in Fig. 0. Notice that these fluxes are extremes and do not represent 
typical errors. Since the direct observation of cosmic-ray composition is limited to energies less 
than 10 6 GeV/nucleus, and there may be a possible steepening of the proton flux above 40 TeV 
(Asakimori et al. 1993), flux of diffuse gamma-rays above about 10 4 GeV suffers large uncertainty. 



2.4. Effect of heavier nuclei 

The effect of heavier nuclei than proton in cosmic-ray and the interstellar matter was estimated 
by different authors. It has been treated as a constant multiplification factor independent of 
energy. This is valid for the calculation of diffuse gamma-rays in the GeV region, since the 
proton/helium ratio in cosmic-rays seems to be constant in the 10~100 GeV range. However, 
recent observations show different spectral indices between protons and other heavier nuclei (see 
Biermann, Gaisser and Stanev 1994; Wiebel-Sooth, Bierman and Meyer 1995 for summary). This 
difference may be explained by a model in which protons and heavier nuclei come from different 
kinds of sources (Biermann, Gaisser and Stanev 1994). Therefore we derive an energy dependence 
of this nuclear enhancement factor, e M , assuming E~ 2 ' 75 spectrum for protons and j? -2 - 63 spectra 
for heavier nuclei (Wiebel-Sooth, Bierman and Meyer 1995). Following Gaisser and Shaefer (1992), 
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e M can be written as 



M 




&ip , @ia I flat 

Wip 1 — J "V 

Cpp Cpp \ re H/ISM 



(4) 



dropping a factor which is related to propagation (only necessary for the case of antiprotons 
discussed in Gaisser and Shaefer (1992), where (rii/n p )cR is the ratio of numbers of nuclei of type 
i to protons (the subscript CR indicates the quantity in cosmic-rays), Wi p {wi a ) the total number of 
wounded nucleons in a collision between a nucleus i and a proton (a), (Ji p {ai a ) the total inelastic 
cross section for a collision between a nucleus i and a proton (a), and (n a /nu)isM is the ratio 
of helium to hydrogen in interstellar matter (ISM). The ISM is assumed to be a mixture of 93% 
hydrogen and 7% helium after Garcia- Munoz et al. (1987). We assume 



?(T P ) 



71, 



^i(GeV) 



CR 



n. 



CR 



( T P ) 

\ 100 GeV J 



0.12 



(i = p or T p < 100 GeV) 
(otherwise) . 



(5) 



where (nj/n p (GeV))cR is the cosmic-ray abundance in GeV range and is tabulated in Gaisser and 
Shaefer (1992). The result is plotted in Fig. || which is to be compared with previous estimates: 
1.5 by Cavallo and Gould (1971), 1.6 ± 0.1 by Stephens and Badhwar (1981), 1.45 by Dermer 
(1986a). e M takes a value of 1.52 for T p < 100 GeV, as given in Gaisser and Shaefer (1992). Notice 
that in the above expression the relativistic rise of interaction cross sections is cancelled out. The 
differences may be partly attributed to a larger multiplicity enhancement factor used in Stephens 
and Badhwar (1981) and a lower helium/proton ratio used in Dermer (1986a). 



2.5. Monte Carlo calculation 

Monte Carlo simulators can produce gamma-rays as final products instead of 7r°'s, including 
the decay kinematics. Also other secondaries, such as neutral kaons which may yield gamma-rays, 
are included, although their contribution is minor. Thus we generated histograms for gamma-rays, 
not for 7r°'s, for various proton energies with logarithmic binning (1 decade = 10 bins). In this 
case the equation ([!]) can be written in a summation form of 

q{E ld ) ■ AEj = 4im H ^ iST Pii j p (T p ^{Ca^(T p ^)f{T Pji \ E 7tj ), (6) 

i 

or, using AE 7>j = E yJ A and AT Pji = T P;i A where A = 10 a05 - 10~ a05 , 

l( E i,j) = ^nH^2jp(T Pti )(C(T w (T Pt i))f(T P) i) Tpj/E^j (7) 

i 

with f(T Pt i;E 7t j) the normalized value of the histograms at Eyj (J2i /(^p,ij E~j) = 2 here to 
incorporate the fact that 27's are produced by a ir°). In the range of energies 3 GeV < T p < 12.5 
GeV no complete experimental data are available. Therefore we decided to use Hadrin for T p < 8 
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GeV, the scaling model for T p > 3 GeV, and Fritiof for T p > 10 GeV: we joined models with a 
linear connection in the overlapping energy regions (3-8 GeV) to obtain a smooth curve. Above 
10 GeV results from two models will be compared. For Monte Carlo programs 80,000^1,000 
events were generated for incident proton energies of T p = 10 -0 ' 5 , 10~ ' 4 , 10 -0 ' 3 , . . . , 10 69 , 10 
GeV, respectively, depending on energy in order to have enough statistical accuracy while saving 
computing time. 

3. Results 

Fig. |9] shows our results on the emissivity of gamma-rays from the interaction of cosmic-rays 
with unit density of atomic hydrogen for two interaction models for T p > 10 GeV, the scaling 
model and Fritiof. The statistical errors due to Monte Carlo calculations are (typically): 4, 0.6, 
0.8, 2, and 3% at E y = 0.01, 0.1, 1, 10, and 100 GeV, respectively. For the cosmic-ray proton flux, 
we assumed three cases: "median" , "maximum" , and "minimum" . (For Fritiof we only show the 
result with the "median" flux.) The "maximal" ("minimal") proton flux gives about 50% higher 
(lower) gamma-ray flux at around 0.1 GeV and about 20% at around 10 GeV compared with that 
with the "median" flux. Remembering that these are extremes, we may say that the standard 
uncertainty of the calculated flux derived from the uncertainty of cosmic-ray proton flux is around 
20 to 30%. The calculated emissivity is 34% higher at 10 GeV when we use Fritiof compared 
with the scaling model: this may be a measure of systematic uncertainty related to the interaction 
model of the present predicition. 

Table @ summarizes the calculated gamma-ray emissivities for three cases of the cosmic-ray 
proton flux in differential and integral form with the scaling model. 
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Table 1. Gamma-ray emissivity due to secondary particle production in collisions of cosmic-rays 

with interstellar matter. 



Energy Differential rate in Integral rate in 

in photon s _1 GcV _1 n IJ 1 photon s -1 ™^ 1 

GeV minimum median maximum minimum median maximum 



1.000E-02 


5 


.872E- 


■26 


9 


.960E- 


■26 


1. 


380E- 


■25 


9. 


.335E-26 


1 


.646E-25 


2 


.385E- 


■25 


1.585E-02 


1 


.133E- 


■25 


2 


.019E- 


-25 


2. 


924E- 


-25 


9, 


.298E-26 


1 


.639E-25 


2 


376E- 


■25 


2.512E-02 


1 


.894E- 


■25 


3 


.540E- 


■25 


5, 


358E- 


-25 


9 


.188E-26 


1 


.619E-25 


2 


.347E- 


■25 


3.981E-02 


2 


.478E- 


■25 


4 


.726E- 


-25 


7. 


298E- 


-25 


8, 


.915E-26 


1 


.568E-25 


2. 


.268E- 


■25 


6.310E-02 


2 


.721E- 


■25 


5 


.223E- 


■25 


8, 


124E- 


■25 


8. 


^384E-26 


1 


.466E-25 


2 


112E- 


■25 


1.000E-01 


2 


.581E- 


■25 


4 


.937E- 


■25 


7, 


652E- 


-25 


7. 


.495E-26 


1 


.296E-25 


1 


.847E- 


■25 


1.585E-01 


2 


.088E- 


■25 


3 


.937E- 


-25 


6. 


006E- 


■25 


6. 


206E-26 


1 


.050E-25 


1 


.467E- 


-25 


2.512E-01 


1 


.333E- 


■25 


2 


.413E- 


-25 


3, 


.549E- 


-25 


4, 


.651E-26 


7 


.594E-26 


1 


.027E- 


■25 


3.981E-01 


6 


.908E 


■26 


1 


.175E- 


■25 


1. 


630E- 


-25 


3. 


164E-26 


4 


.940E-26 


6 


.416E- 


■26 


6.310E-01 


3 


.206E- 


■26 


5 


.116E- 


-26 


6. 


684E- 


-26 


1 


.977E-26 


2 


.949E-26 


3 


690E- 


■26 


1.000E+00 


1 


301E- 


■26 


1 


961E- 


-26 


2. 


447E- 


-26 


1. 


.131E-26 


1 


.618E-26 


1 


.971E- 


■26 


1.585E+00 


4 


.765E- 


-27 


6 


.861E- 


-27 


8, 


.314E- 


-27 


5, 


.986E-27 


8 


.229E-27 


9 


.847E- 


-27 


2.512E+00 


1 


.546E- 


■27 


2 


125E- 


-27 


2, 


.524E- 


-27 


2. 


.982E-27 


3 


.946E-27 


4 


.681E- 


■27 


3.981E+00 


4 


.812E- 


-28 


G 


.348E- 


-28 


7. 


.494E- 


-28 


1. 


.455E-27 


1 


.866E-27 


2 


.216E- 


-27 


6.310E+00 


1 


.497E- 


-28 


1 


.910E- 


-28 


2. 


266E- 


-28 


7. 


.027E-28 


8 


.805E-28 


1 


052E- 


■27 


1.000E+01 


4 


.579E- 


■29 


5 


.698E 


■29 


6. 


.822E- 


-29 


3. 


.324E-28 


4 


.108E-28 


4 


934E- 


■28 


1.585E+01 


1 


.361E- 


■29 


1 


.673E- 


■29 


2. 


016E- 


-29 


1 


539E-28 


1 


.894E-28 


2 


.280E- 


■28 


2.512E+01 


3 


.952E- 


■30 


4 


.854E- 


-30 


5, 


.852E- 


-30 


7. 


028E-29 


8 


.672E-29 


1 


.041E- 


-28 


3.981E+01 


1 


.133E- 


-30 


1 


.398E- 


30 


1. 


678E- 


-30 


3 


192E-29 


3 


.955E-29 


4 


.731E- 


■29 
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At this stage it may be worth showing the contributions of cosmic-ray protons of various 



energies to gamma-rays of some specific energies in Fig. 10, From this figure one can see most 
diffuse gamma-rays below about 300 MeV of nuclear origin come from protons with kinetic 
energies of a few to several GeV, where no experimental results for secondary pion production are 
available and where there is a large uncertainty in the cosmic-ray proton flux. 

The results obtained here are compared with previous calculations (Cavallo and Gould 1971, 
Stecker 1979, Stephens and Badhwar 1981, Dermer 1986a) in Figs. |ll] and |l2[ where in the latter 
they are plotted in integrated form: q(> E^) = JJ° q(Ei f )dE! y . The main difference between 
authors may come from the assumed demodulated spectrum of the cosmic-ray protons. (See Fig. 1 
of Dermer (1986a) for comparison of assumed spectra.) 

We can extend our calculation to higher energy gamma-rays and the results are presented 
in Fig. [l3| with previous calculations (Stecker 1979, Berezinsky et al. 1993, Chardonnet et al. 
1995). For this plot we show calculation only with the scaling model: Fritiof causes fatal errors 
at ^ 10 TeV. Since the mean interaction length of gamma-rays for photon-photon collisions with 
the cosmic microwave background radiation becomes comparable to the distance to the Galactic 
center from Earth above several 100 TeV (Protheroe 1986), gamma-ray flux at higher energies 
than this are greatly reduced. Thus we stop our calculation at 100 TeV. The proton fluxes used 
in other calculations are: 2.35£"~ 2 ' 67 (Stecker) and 1.75-E^ 2 ' 73 (Chardonnet et al.), where E p is in 
GeV and the unit is cm" 2 s _1 sr _1 GeV _1 . Berezinsky et al. assumed the total cosmic-ray spectrum 
to be 1.59£p~ 2 ' 73 cm~ 2 s~ 1 sr~ 1 (GeV/nucleon) _1 . Also shown is the effect of energy-dependent 
nuclear enhancement factor (e M ), which raise the expected flux by about 10% at 1 TeV. The flux 
given in Berezinsky et al. is smaller than others, partly due to the assumed cosmic-ray spectrum. 

In Tabled the integral gamma-ray emissivities calculated by different authors are summarized. 



4. Discussion 

Now we compare our results with the EGRET data (Hunter et al. 1996). 

Fig. [TJI shows the Galactic diffuse gamma-ray flux in a differential form, multiplied by E^ 
to make the differences more visible. The dotted line shows the expectation as the sum of three 
components of diffuse gamma-rays: cosmic-ray nuclear interactions (which is the main concern 
of this paper and is dominant above about 200 MeV), electron bremsstrahlung, and inverse 
Compton. The components were calculated by the model developed by Hunter et al. (1996), which 
incorporates the detailed distribution of interstellar matter, and are averaged over the Galactic 
center region, 300° < £ < 60°, |6| < 10°. In their model they used the gamma-ray emissivity from 
cosmic-ray interactions given by Stecker (1988). Solid and dashed lines are the expected fluxes 
in which the gamma-ray emissivity is replaced with the ones computed in this work: solid line 
assumes the "median" proton flux and dashed lines assume the "minimal" and "maximal" fluxes. 
None of the models mentioned above seems to yield enough gamma-rays above 1 GeV: thus we 
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Table 1 — Continued 



Energy Differential rate in Integral rate in 

in photon s _1 GcV^ 1 n^ 1 photon s^ 1 n }i 1 

GeV minimum median maximum minimum median maximum 



G 


.310E+01 


3 


.238E- 


■31 


4 


.015E- 


■31 


4, 


799E- 


■31 


1. 


.451E- 


-29 


1 


.804E-29 


2 


151E- 


■29 


1. 


.000E+02 


9 


.277E- 


32 


1 


.153E- 


■31 


1. 


375E- 


■31 


6 


.625E- 


■30 


8. 


.248E-30 


9 


.820E- 


■30 


1. 


.585E+02 


2 


.669E- 


■32 


3 


.324E- 


-32 


3, 


.956E- 


-32 


3. 


039E- 


■30 


3 


.788E-30 


4 


505E- 


■30 


2 


.512E+02 


7. 


.721E- 


■33 


9 


.626E- 


■33 


1. 


145E- 


32 


1 


.403E- 


-30 


1 


.750E-30 


2 


.079E- 


30 


3 


.981E+02 


2 


.248E- 


33 


2 


.805E- 


■33 


3, 


333E- 


■33 


6 


.514E- 


■31 


8. 


129E-31 


9 


.656E- 


■31 


6 


.310E+02 


6 


.590E- 


■34 


8 


.224E- 


-34 


9, 


.768E- 


-34 


3. 


.042E- 


■31 


3 


.797E-31 


4 


509E- 


■31 


1 


.000E+03 


1 


.941E- 


■34 


2 


.423E- 


-34 


2, 


.878E- 


-34 


1 


.427E- 


■31 


1 


.782E-31 


2 


115E- 


■31 


1 


.585E+03 


5 


.747E- 


■35 


7 


.175E- 


-35 


8, 


.518E- 


-35 


6 


724E- 


■32 


8. 


.395E-32 


9 


.967E- 


■32 


2 


.512E+03 


1 


.708E 


■35 


2 


.133E- 


-35 


2. 


532E- 


-35 


3 


180E- 


-32 


3 


.971E-32 


4 


.714E- 


■32 


3 


.981E+03 


5 


.098E- 


■36 


6 


.365E- 


-36 


7, 


556E- 


-36 


1 


.510E- 


-32 


1 


.885E-32 


2 


238E- 


32 


6 


.310E+03 


1 


.527E- 


■36 


1 


.906E- 


-36 


2, 


263E- 


-36 


7. 


188E- 


■33 


8. 


.976E-33 


1 


.065E- 


■32 


1 


.000E+04 


4 


.587E- 


■37 


5 


.727E- 


■37 


6, 


798E- 


37 


3 


.432E- 


■33 


4 


.285E-33 


5 


.086E- 


33 


1 


.585E+04 


1 


.383E- 


■37 


1 


.726E- 


■37 


2. 


049E- 


-37 


1 


642E- 


■33 


2 


.050E-33 


2 


434E- 


33 


2 


.512E+04 


4 


.181E- 


■38 


5 


.221E- 


■38 


6. 


197E- 


■38 


7. 


863E- 


■34 


9 


.819E-34 


1 


165E- 


33 


3 


.981E+04 


1 


.269E- 


■38 


1 


.585E- 


■38 


1. 


.881E- 


-38 


3. 


.759E- 


-34 


4 


.694E-34 


5 


.572E- 


■34 


6 


.310E+04 


3 


.856E- 


■39 


4 


.816E- 


-39 


5. 


716E- 


-39 


1 


.784E- 


-34 


2 


.227E-34 


2 


644E- 


■34 


1 


.000E+05 


1 


153E- 


■39 


1 


.440E- 


-39 


1. 


.709E- 


39 


8. 


340E- 


-35 


1 


.042E-34 


1 


.236E- 


-34 



Note. — Figures to be read as 5.782E-26 = 5.782 x 10~ 26 , etc. "Minimum", "median" and 
"maximum" denote the models of the cosmic-ray spectrum (see text). The scaling model is used 
for T p > 10 GeV. 
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may conclude that even the detailed model studied here cannot explain the excess of Galactic 
diffuse gamma-ray flux above 1 GeV. 



As one of the possible solutions to the problem, we tried to fit the cosmic-ray proton spectrum 

p( e p) = aE p' 



from the EGRET data with a single power-law spectrum in total energy, J p (E p ) = aE b , and 



have computed the x 2 ' s i n the a-b plane and searched the minimum. 

The results are a = 1.481q2o cm _2 s _1 sr _1 GeV _1 and b = 2.451qo4 when we use the scaling 
model for T p > 10 GeV (a = 1.67t°;!o cm-W^GeY- 1 and b = 2.51 ± 0.04 with Fritiof for 



T p > 10 GeV). The expected flux calculated from this best fit is plotted in Fig. |l4j as a dotdashed 
line. This spectrum gives a better fit above 1 GeV, but predicts higher gamma-ray flux than 
observed around several 100 MeV. 

However, this flatter proton spectrum is not consistent with the direct observations above 
10 GeV (see Fig. ^). Therefore one may be tempted to suggest a possibility that the cosmic-ray 
spectrum is flatter in the Galactic center region than the one in the local (solar neighborhood) 
region. It is interesting to note that the power-law index obtained here is similar to that of 
the cosmic-ray source spectrum predicted by the cosmic-ray reacceleration theory (Ptuskin 
1995 and references therein): in this theory the source spectrum is proportional to i?~ 7s , where 
7 S = 2.3 ~ 2.4 (where R is the rigidity) and the escape length to i? -1 / 3 to match the observed 
spectrum at all energies. 

In any case, more detailed analysis is necessary before going on, since in the above fitting 
the bremsstrahlung and inverse Compton components are fixed: these components will change 
following the electron spectrum which is related to the proton spectrum through the cosmic-ray 
proton to electron ratio, which is usually believed to be universal. 

On the other hand, COS B data indicate that the diffuse gamma-ray spectrum in the inner 
Galaxy (310° < I < 50°) is consistent with the cosmic-ray proton spectrum of about cx E~ 2 ' 7 
(Bloemen 1987), which does not seem to be consistent with the EGRET result (Hunter et al. 
1996). 

Further observations, including analysis of the EGRET data at high latitudes (Sreekumar et 
al. 1996), are nessesary to resolve the spectral variation problem. 



5. Conclusion 

The gamma-ray spectrum from the decay of secondary neutral pions produced by interactions 
of cosmic-ray protons with interstellar gas has been calculated utilizing modern event generators 
which simulate high-energy p-p collisions by Monte Carlo methods. The result is not inconsistent 
with previous calculations: the observed excess of the Galactic diffuse gamma-rays above about 
1 GeV can not be explained by revising the models of cosmic-ray interaction with interstellar 
matter. This might suggest flatter cosmic-ray proton spectra in the central region of the Galaxy 
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than in the solar neighorhood, but requires further study. 

I gratefully acknowledge the EGRET team for the kind hospitality during my stay at 
NASA/GSFC, especially Dr. Stan Hunter for discussions and providing me the EGRET diffuse 
gamma-ray fluxes prior to publication, and Dr. Robert Hartman for his continuous encouragement 
and support. I also thank Drs. Jonathan Ormes, Charles Dermer, and Felix Aharonian for valuable 
discussions and suggestions. 
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Table 2. Integrated Galactic gamma-ray emissivity from cosmic-ray interactions with interstellar 

matter. 



Reference 


Model 






q(> E,r 







0.1 GeV 


1 GeV 


1 TeV 


Cavallo and Gould 1971 


Isobar + Phase space 


2.4 


1.8 


0.13 




Stephens and Badhwar 1981 


Scaling 


1.92-2.34 


1.37-1.63 


0.115-0.116 




Dernier 1986a 


Isobar + Scaling 


2.02 


1.53 


0.159 




Stecker 1988 b 


Isobar + Scaling 


1.97 


1.59 


0.209 




Berezinsky et al. 1993 


SIBYLL 








0.76 x 10~ 6 


Chardonnet et al. 1995 


Pythia 








1.5 x 10~ 6 


This work c 


Hadrin + Scaling 


1.65 


1.30 


0.162 


1.78 x 10~ fi 






(0.94-2.39) 


(0.75-1.85) 


(0.113-0.197) 


((1.43-2. 12)xl0~ 6 ) 



a Units: 10 25 s 1 per hydrogen atom. 

b Computed using the parameterization given in Bertsch et al. (1993). 

c The "median" proton flux is used. (Emissivities with the "minimal" and "maximal" proton fluxes are shown in 
parentheses.) 
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Fig. 1. — Comparison of the model predictions with experimental data of the secondary 7r° spectrum 
in the CMS, corresponding to LS proton kinetic energy T p = 0.97 GeV. Solid histogram: Hadrin, 
dashed line: isobar model by Stecker (1970), dotted line: scaling model by Stephens and Badhwar 
(1981), and asterisks: data of Bugg et al. (1964). 

Fig. 2. — Same as in Fig. [I], but with the secondary ir° spectrum in the laboratory system (LS). 

Fig. 3. — Comparison of the model prediction with experimental data of the rapidity distribution 
of secondary gamma-rays at LS proton momentum of 12.4, 205 and 300 GeV/c. Solid (dotted) 
histograms: Fritiof (Pythia), dashed curves: combination of the scaling model by Stephens & 
Badhwar (1981) and Decay. Calculations are normalized to the inclusive cross section given by 
Dermer (1986b). Data of Jaeger et al. (1975a), Jaeger et al. (1975b) and Sheng et al. (1975). 

Fig. 4. — Comparison of the invariant cross section for the production of neutral pions as a function 
of pion momentum in CMS at different emission angles for various colliding beam energies. Data 
points and scaling calculation (solid lines) are taken from Stephens and Badhwar (1981). Solid 
(dotted) histograms are converted from Fritiof (Pythia) results using the Sternheimer relation 
(see text). 

Fig. 5. — Inclusive cross sections for the production of neutral pions in proton-proton collisions as a 
function of the incident proton momentum. The data points are from the compilation of Appendix 
A of Dermer (1986b) and the dashed lines are the fits by Dermer (1986b). The results of Fritiof 
(Pythia) are plotted in dotted (dotdashed) lines and a fit used in this work is shown by a solid 
line. 

Fig. 6. — Compilation of the demodulated cosmic-ray proton spectra calculated by different authors 
and the data on the cosmic-ray proton flux at high energies. Lines are the demodulated spectra 
(Gloeckler and Jokipii 1967, Stephens and Badhwar 1981, Simpson 1983, Ormes and Protheroe 
1983, Ip and Axford 1985, Nagashima et al. 1989, Webber and Potgieter 1989, Seo et al. 1991) and 
symbols are the data at high energies (pluses, open circles, closed circles, squares, crosses, inverse 
triangles and triangles are from Ryan et al. 1972, Smith et al. 1973, Webber et al. 1987, Kawamura 
et al. 1989, Seo et al. 1991, Asakimori et al. 1993 and Ivanenko et al. 1993, respectively.) Also 
shown by thick solid line is the "median" flux used in this work. See also Fig. 0. (The upper and 
lower dashed lines correspond to the curves BS and Mu of Stephens and Badhwar (1981).) 

Fig. 7. — The demodulated cosmic-ray proton spectra used in this work: "median", "maximum", 
"minimum" and "best fit", where the "maximum" ("minimum") is taken as the upper (lower) 
envelopes of different calculations and data shown in Fig. || "Best fit" is derived from the EGRET 
data (see text). Also shown are the spectra used in previous calculation (Stephens and Badhwar 
1981, Dermer 1986a) of the gamma-ray emissivity. (The upper and lower dashed lines correspond 
to the curves BS and Mu of Stephens and Badhwar (1981).) 

Fig. 8. — The nuclear enhancement factor, e M , as a function of the incident cosmic-ray proton 
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energy. 

Fig. 9. — The gamma-ray production source function per unit density of atomic hydrogen from the 
interaction of the cosmic-ray with the interstellar medium for different models of high-energy p-p 
collisions. Solid and dotted lines are the results with the scaling model and Fritiof, respectively, 
for T p > 10 GeV assuming the "median" cosmic-ray proton flux. Both of them used the same model 
for T p < 10 GeV (see text). Also shown are the results assuming the "minimal" and "maximal" 
proton fluxes with the scaling model. 

Fig. 10. — The contributions of cosmic-ray protons of various energies to gamma-rays of some 
specific energies. 

Fig. 11. — Comparison of the gamma-ray production source functions per unit density of atomic 
hydrogen from the interaction of the cosmic-ray with the interstellar medium by different authors. 
Triple-dot-dashed, dashed, dotdashed lines are differential emissivities multiplied by E 2 given in 
Cavallo and Gould (1971), Stephens and Badhwar (1981), Dermer (1986a) and Stecker (1988), 
respectively (the dotted curve is drawn by the parameterization given in Bersch et al. (1993), 
which fits the results of Stecker (1988)), and the solid line is from this work with Pythia and the 
"median" proton flux. 

Fig. 12. — Comparison of the integrated gamma-ray production source functions per unit density 
of atomic hydrogen from the interaction of the cosmic-ray with the interstellar medium by different 
authors. Notations are the same as in Fig. [| 

Fig. 13. — The gamma-ray production source functions per unit density of atomic hydrogen from 
the interaction of the cosmic-ray with the interstellar medium for different models of high-energy 
p-p collisions at high energies. Solid lines are the results with the scaling model and the "median" 
proton flux, dashed lines are similar but with a fixed nuclear enhancement factor (e M = 1.52). Also 
shown are previous calculations (Stecker 1979, Berezinsky et al. 1993, Chardonnet et al. 1995). 

Fig. 14. — The diffuse gamma-ray differential spectrum, multiplied by E 2 , of the Galactic center 
region. Data are from EGRET observations (Hunter et al. 1996) averaged over 300° < £ < 60°, 
\b\ < 10°. The model predictions based on gamma-ray emissivities from cosmic-ray nuclear 
interactions with different assumptions are also shown. Solid line: the scaling model and the 
"median" proton flux, dashed lines: the scaling model and the "maximum" and "minimum" proton 
flux, and dotted line: from Hunter et al. (1996) based on the result of Stecker (1988). Also shown 
by dotdashed lines are spectra obtained with the scaling model and the "best fit" proton flux 
(oc Ep 2 ' 45 ). See text for detail. 
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